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1.  INTRODUCTION:  Research  Objectives  and  Summary  of  the  Results 

1.1  Aims  of  the  Research 

The  general  objective  of  this  research  was  to  determine  by  theoretical  simulations 
dynamical  properties  of  solid  hydrogen  clusters,  doped  with  certain  atoms,  such  as  boron, 
magnesium  and  lithium.  Doped  solid  hydrogen  systems  are  of  great  interest  as  potentially 
advanced  propellants  of  high  Specific  Impulse.  (Ref.  1)  As  such,  knowledge  of  the  prop¬ 
erties  of  such  systems  is  of  major  importance  for  the  Air  Force  project  on  High  Energy 
Density  Materials  (HEDM).  Also,  doped  solid  hydrogen  and  doped  hydrogen  clusters  are 
of  basic  scientific  interest,  since  such  systems  are  expected  to  exhibit  strong  quantum 
mechanical  effects.  The  properties  targeted  in  the  proposal  involved  time-dependent  phe¬ 
nomena  and  processes.  An  objective  of  central  importance  of  the  research  was  the  devel¬ 
opment  of  methods  and  algorithms  for  the  simulation  of  time-dependent  behavior  in  low- 
temperature  systems  exhibiting  strong  effects.  The  challenge  here  is  to  develop  a  novel 
computational  method  for  solving  the  quantum-mechanical  time-dependent  Schrodinger 
equation  for  many-atom  systems.  This  goal,  described  as  Task  7  of  the  original  Program 
Plan,  was  essential  for  constructing  a  tool  of  simulation  for  the  purposes  of  this  proposal, 
but  clearly  if  attained  it  provides  a  powerful  means  also  for  future  studies  of  many  other 
properties  of  low-temperature,  highly  quantum-mechanical  systems.  One  of  the  applica¬ 
tions  targeted  in  this  proposal  was  the  description  of  the  dynamics  following  excitation  of 
a  metal  atom  in  a  solid  hydrogen  environment,  or  in  hydrogen  clusters.  This  is  important 
for  the  spectroscopic  characterization  of  the  system.  This  is  the  content  of  Tasks  6  and  8 
of  the  original  proposal.  Another  important  time-dependent  property  is  the  dynamical  sta¬ 
bility  of  the  system.  One  of  the  objectives  of  the  proposal  (Task  5  of  the  revised  Work 
Plan)  was  to  calculate  collisions  of  doped  hydrogen  clusters  with  atoms  and  molecules, 
and  determine  their  survival  probability.  This  information  is  important  for  some  of  the 
methods  for  preparing  doped  hydrogen  systems.  To  fully  characterize  doped  hydrogen 
systems,  structural  and  energetic  properties  must  be  described,  in  addition  to  the  time- 
dependent  behavior.  One  of  the  key  objectives  of  the  proposal  was  to  determine  by  quan¬ 
tum-mechanical  simulations  the  structural  and  energetic  properties  of  various  doped 
hydrogen  clusters  (Task  1  and  Task  4  of  the  original  workplan).  Another  objective  was  to 
determine  the  stability  of  doped  hydrogen  systems  against  diffusion  (Task  3  of  the  Work 
Plan).  In  summary,  the  main  aims  of  this  proposal  were  to  develop  simulation  methods 
for  many-atom  quantum  systems,  and  to  carry  out  calculations  for  several  time-dependent 
properties,  as  well  as  for  related  structural  and  energetic  properties  of  doped  hydrogen 
systems. 

1.2  Brief  Summary  of  the  Research  Results 

In  the  course  of  the  research  project  described  in  this  report,  the  following  main 
results  were  obtained: 

(1)  A  novel  method  was  developed  for  time-dependent  quantum  simulations  of 
large  systems.  The  method  is,  to  the  author’s  knowledge,  a  unique  tool  to  date  for  quan¬ 
tum  simulations  in  time  of  processes  in  many-atom  systems.  Applications  for  realistic 
systems  having  up  to  -  50  atoms  are  at  hand,  extensions  to  much  larger  systems  are  under 
way.  The  method  has  a  wide  range  of  potential  applications  for  HEDM  systems,  e.g. 
doped  solid  hydrogen,  for  cryogenic  systems  in  general,  and  also  for  other  molecular  sys¬ 
tems.  Preliminary  applications  are  described  in  the  report.  (Refs.  2-4)  This  development 
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successfully  accomplishes  Task  7  of  the  original  Program  Plan.The  author  regards  this  as 
the  most  important  result  of  this  project. 

(2)  Another  time-dependent  simulation  method,  based  on  the  Time-Dependent  Self- 
Consistent  Field  (TDSCF)  approximation,  was  also  developed  and  tested  for  quantum 
clusters.  Not  as  efficient  and  widely  applicable  in  general  as  method  (1)  above,  the 
TDSCF-based  algorithm  was  found  to  have  a  useful  "niche"  of  applications  for  cryogenic 
materials  where  quantum  effects  are  moderate,  such  as  solid  Ne,  or  solid  H2  at  high  pres¬ 
sures.  This  result  also  contributes  to  the  accomplishment  of  Task  7  of  the  Program  Plan. 

(3)  The  structural  and  energetic  properties  of  B(H2)„  clusters,  n  <  8,  were  studied  by  Dif¬ 
fusion  Quantum  Monte  Carlo  (DQMC)  simulations  (cooperation  with  Prof.  M.H.  Alexan¬ 
der).  For  the  smaller  clusters,  the  structures  are  planar,  with  the  H2  molecules  located  nor¬ 
mal  to  the  p  orbital  of  the  boron  atom.  The  B  atom  is  solvated  inside  a  ring  of  H2 
molecules,  for  n  >  6;  that  is,  the  boron  is  encircled  by  the  hydrogen,  and  stabilized  in  this 
way.  The  binding  energy  is  much  higher  than  for  neat  (H2)„  clusters.  In  the  solid  one 
therefore  expects  a  strong  cluster  involving  the  B  atom  and  at  least  6  neighboring  H2 
molecules.  This  should  be  stable  up  to  the  melting  point  of  solid  hydrogen.  Since  the  B 
atom  is  strongly  bound  to  neighboring  H2  molecules,  the  diffusivity  of  B  in  solid  H2  is 
predicted  to  be  extremely  low.  This  (together  with  corresponding  results  for  Mg(tf2)„ 
clusters)  accomplishes  Tasks  1  and  4  of  the  updated  Program  Plan. 

(4)  The  predicted  "solvation"  of  B  atoms  in  H2,  and  the  strong  clustering  of  atoms  with 
several  neighboring  H2  molecules  with  resulting  low  diffusivity  even  for  temperatures 
well  above  4K,  suggest  that  an  appreciable  concentration  of  B  atoms  (several  percent)  in 
solid  hydrogen  should  be  fairly  stable.  By  these  results,  boron  in  hydrogen  appears  a  very 
promising  HEDM  material,  far  superior  to  any  other  system  studied  in  this  project. 

(5)  Structural  and  energetic  properties  of  Mg(H2)„  clusters,  with  n  <  18,  were  studied  by 
DQMC.  In  the  ground  state,  which  the  DQMC  computes,  the  magnesium  is  solvated 
inside  a  shell  of  H2  molecules.  The  shell  is  not  isotropic  for  any  cluster  size,  and  can  be 
roughly  described  as  an  anisotropic  quantum  liquid.  However,  the  Mg  atom  cannot  be 
found  at  or  beyond  the  H2  shell  with  any  significant  probability,  so  it  is  "shielded"  inside 
the  H2  shell.  This  (together  with  the  above  results  for  boron,  and  with  (6)  below)  accom¬ 
plishes  Tasks  1  and  4. 

(6)  Using  Feynman  Path  Simulations,  structural  properties  of  Mg(H2)„  at  thermal  equilib¬ 
rium  were  calculated  for  different  temperatures.  It  was  found  that  the  Mg  remains  sol¬ 
vated  inside  the  shell  of  H2,  for  12  <  n,  up  to  the  evaporation  point  of  the  cluster,  around 
4.2K.  On  this  basis,  diffusion  of  Mg  in  solid  hydrogen  is  expected  to  be  very  low  under 
4K.  This,  together  with  previous  results  above,  accomplishes  Tasks  1, 4. 

(7)  Mg  in  hydrogen  has  properties  useful  for  HEDM  by  our  results,  though  not  as  advan¬ 
tageous  as  those  of  B  in  H2.  The  low  diffusivity  below  4K,  due  to  the  solvation  effect, 
suggests  that  Mg  in  appreciable  concentrations  in  hydrogen  can  remain  stable  at  low  tem¬ 
peratures.  Mg  seems  much  superior  in  this  respect  to  Lithium  in  hydrogen. 

(8)  The  ground-state  vibrational  wavefunction  of  Mg(H2)„  is  highly  delocalized.  The 
Mg-H2  distance  distribution  covers  relative  distances  from  -2. 5 A  to  ~7.  oA.  In  the 
larger  clusters,  n  >  10,  the  computed  wavefunction  suggests  "radial"  vibrations  of  the  H2 
with  respect  to  the  Mg,  and  collective  vibrations  of  the  hydrogens  on  a  near-spherical 
shell.  The  characterization  of  the  vibrational  wavefunction  of  Mg(H2)„  accomplishes 
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(together  with  (9)  below),  Task  2. 

(9)  The  ground-state  vibrational  wavefunction  of  B(H2)„  is  also  very  delocalized, 
although  the  relative  B  -  H2  distance  distribution  is  less  spread  out  than  the  Mg  -  H2  one 
(the  range  of  distances  is  3 A  -5 A,  roughly).  The  mean  structure,  say  of  B(H2)6,  is  pla¬ 
nar,  but  the  out  of  plane  zero-point  structural  fluctuations  are  very  large.  This  characteri¬ 
zation  of  the  vibrational  motion  accomplishes  part  of  Task  2  of  the  Program  Plan. 

(10)  Using  time-dependent  quantum  mechanical  methods,  collisions  of  quantum  clusters 
with  atoms  or  with  molecules  were  studied.  A  surprising  and  potentially  important  find¬ 
ing  was  that  survival  probabilities  of  the  quantum  clusters  in  such  collisions  are  quite 
high,  despite  the  low  binding  energies.  For  instance,  at  a  collision  energy  of  150K  the 
survival  probability  of  a  B  -  H2  cluster  hit  by  an  H2  molecule  is  about  50%  !  These 
results  on  survival  of  quantum  clusters  in  collisional  impact  fulfill  Task  5  of  the  updated 
Workplan.  The  results  point  to  possible  strategies  for  the  preparation  of  HEDM  materials 
(such  as  boron  in  solid  H2),  based  on  deposition  of  quantum  clusters. 

(11)  Using  a  time-dependent  quantum  algorithm,  we  calculated  the  dynamics  following 
electronic  excitation  of  small  quantum  clusters  such  as  Li(H2)2.  A  long  lifetime,  of  the 
order  of  104  femtoseconds,  was  found  for  the  excited  state.  The  channel  producing 
Li(P)  +  (H2)2  was  found  to  dominate  over  the  channel  resulting  in  Li(P)  -  H2  +  H2.  The 
results  have  applications  to  characterization  of  quantum  clusters  (of  HEDM  interest)  by 
electronic  spectroscopy,  including  by  ultrafast  spectroscopic  techniques.  This  calculation 
attains  the  goal  of  Task  6. 

(12)  The  dynamics  following  electronic  excitation  of  Mg  atoms  in  solid  hydrogen  was 
studied  by  time-dependent  quantum  simulations,  using  the  novel  method  developed  in 
this  project  (item  (1)  above).  The  results  show  coherent  cage  vibrations,  over  a  time- 
duration  of  at  least  0.5  fs.  Photoinduced  mobility  of  Mg  atoms  leaving  the  cage  was 
found.  Spectroscopic  observables  were  computed.  This  calculation  accomplishes  Task  2, 
and  a  large  part  of  Task  8.  It  demonstrates  the  power  of  Quantum  Molecular  Dynamics 
simulations  using  the  new  method. 

(13)  Using  semiclassical  Molecular  Dynamics  methods,  calculations  were  carried  out  on 
the  coupled  electronic-nuclear  motions  for  Fluorine  atoms  in  solid  Kr.  This  cryogenic 
system  is  not  of  direct  HEDM  relevance,  but  is  a  useful  model  for  P-state  open  shell 
atoms  in  solid  H2.  Both  the  dynamics  at  equilibrium,  and  the  dynamics  in  a  photoinduced 
process,  were  simulated.  The  timescales  for  transitions  between  electronic  states,  and  for 
change  of  electronic  polarization  (related  to  electronic  orbital  reorientation)  were  com¬ 
puted.  The  predictions  are  useful  for  spectroscopic  characterization  of  HEDM  systems 
such  as  boron  in  hydrogen.  These  results  contribute  to  the  accomplishments  of  Tasks  2 
and  8. 

(14)  In  cooperation  with  experiments  at  Grenoble,  France,  by  P.  Zeppenfeld,  O.  Wilches 
and  M.  Bienfait,  (Ref.  5)  hydrogen-argon  mixed  monolayers  on  graphite  were  studied. 
This  research,  not  central  to  the  HEDM  issues,  is  nevertheless  important  for  understand¬ 
ing  properties  of  mixtures  in  hydrogen  with  other  atoms.  The  results  show  that  such  mix¬ 
tures  are  inherently  more  disordered  than  corresponding  classical  mixtures  of  particles 
similar  in  sizes  (such  as  Xe+Kr  monolayers).  Also,  the  disordered  structures  seem  to 
depend  on  temperature.  The  results  suggest  a  new  type  of  structural  behavior  of  hydro¬ 
gen- Ar  mixtures  in  two  dimensions,  but  are  still  very  preliminary. 
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With  the  above  results  listed,  all  the  Tasks  set  in  the  original  and  updated  Program 
Plan  have  been  fully  accomplished,  with  the  exception  of  part  of  Task  8.  The  results  add 
significantly  to  our  knowledge  on  doped  hydrogen  clusters  and  solids,  on  cryogenic  sys¬ 
tems  in  general,  and  possible  novel  theoretical  tools  for  the  description  of  processes  in 
such  systems. 

From  a  HEDM  point  of  view,  particularly  important  findings  are  those  which  are 
very  encouraging  as  to  the  suitability  of  boron  in  hydrogen  as  a  potential  HEDM  propel¬ 
lant.  Also  significant  are  the  results  on  the  properties  of  magnesium  in  hydrogen:  They 
point  to  important  advantages  of  this  material  compared  with  other  metal-hydrogen  sys¬ 
tems  hitherto  studied. 


2.  RESEARCH  RESULTS,  DATA  OBTAINED  AND  ACCOMPLISHMENTS 
2.1  New  Method  for  Time-Dependent  Quantum  Simulations  of  Large  Systems 

2.1.1  Introduction:  Classical  trajectory  simulations  have  proved  a  tool  of  pivotal 
importance  for  the  understanding  of  molecular  processes  in  condensed  matter  as  well  as 
in  the  gas  phase.  The  classical  level  of  description  is  certainly  adequate  to  a  large  extent 
for  a  great  variety  of  molecular  phenomena.  This,  together  with  the  computational  sim¬ 
plicity  and  power  of  classical  Molecular  Dynamics  (Ref.  6)  have  made  the  latter  an 
extremely  useful  tool  of  interpretation  and  insight,  and  often  also  of  quantitative  predic¬ 
tion.  However,  in  many  cases  it  is  desirable  to  proceed  beyond  the  classical  level,  to  a 
quantum  description.  Quantum  effects,  e.g.  non-adiabatic  transitions,  zero-point  motion 
effects,  interference  and  tunneling  have  important  experimental  manifestations.  This  is 
relevant  to  many  molecular  systems  at  normal  conditions.  For  systems  such  as  solid 
hydrogen,  quantum  effects  are  very  large,  and  the  ability  to  treat  them  in  simulations  is  of 
critical  importance.  HEDM  systems  such  as  doped  solid  hydrogen  are  clearly  in  this  cate¬ 
gory.  Hence  a  very  central  goal  of  the  research  reported  here  was  to  develop  a  method  for 
time-dependent  quantum-mechanical  simulations  of  systems  having  many  atoms.  This  is 
Task  No.  7  of  the  Program  Plan  of  our  project  -  essentially  to  develop  a  powerful,  practi¬ 
cal  algorithm  for  Quantum  Molecular  Dynamics  applicable  to  large  realistic  systems,  in 
particular  of  the  HEDM  type.  Note  that  despite  impressive  progress,  current  "numerically 
exact"  methods  for  quantum  dynamics  are  confined  to  systems  of  a  few  modes  only  (Ref. 
7).  In  the  author’s  judgement,  the  project  had  important  success  in  producing  a  time- 
dependent  quantum  simulation  method  that  can  treat  large  systems  -  a  major  leap  com¬ 
pared  with  previous  capabilities.  The  method  was  shown  by  tests  to  be  of  good  accuracy, 
and  is  applicable  to  realistic,  complicated  interaction  potentials.  The  novel  method  devel¬ 
oped  is  not  along  the  detailed  lines  envisaged  in  the  original  proposal.  The  Time- 
Dependent  Self-Consistent  Field  (TDSCF)  approximation  (Ref.  8)  and  extensions  thereof 
is  used,  as  a  basis  for  the  Quantum  Molecular  Dynamics  algorithm.  Although  interesting 
results  were  obtained  here  with  the  TDSCF  method  in  this  project,  and  although  in  spe¬ 
cial  "niches"  it  was  shown  to  provide  a  very  useful  approach,  it  cannot  serve  as  the 
"method  of  choice"  for  this  purpose:  The  TDSCF  and  especially  its  extensions,  proved 
computationally  too  costly.  The  new  method,  developed  in  the  present  project,  is  based  on 
the  Classical  Separable  Potential  (CSP)  approximation,  but  includes  also  extensions  of 
high  accuracy,  that  are  nevertheless  computationally  feasible.  The  CSP  method  was 
reported  in  a  series  of  articles  (Refs.  2-4),  though  the  latest  and  most  powerful  extensions 


4 


are  yet  to  be  published.  The  essential  idea  in  the  CSP  is  to  use  a  classical  Molecular 
Dynamics  simulation  for  the  process  of  interest,  to  help  "guide"  the  quantum  simulation, 
and  simplify  the  latter.  This  is  done  by  using  the  classical  MD  to  generate  dynamical 
mean  potentials  for  each  degree  of  freedom.  The  mean  potentials  are  then  used  in  quan¬ 
tum  calculations  to  generate  a  set  of  "orbitals"  for  each  degree  of  freedom  -  orbitals  for 
the  motions  of  atomic  coordinates.  The  time-dependent  multi-dimensional  wavefunction 
of  the  full  system  is  constructed  from  the  "orbitals"  mentioned.  The  SCP  was  tested  and 
developed  first  for  systems  of  moderate  quantum  effects,  e.g.,  cryogenic  heavy  atom  sys¬ 
tems.  However,  the  use  of  classical  MD  only  for  the  limited  purpose  of  constructing 
effective  potentials  also  seems  to  cause  no  difficulties  in  applications  in  solid  hydrogen 
systems.  In  the  following  subsection  the  general  principles  of  the  method  are  described, 
with  applications  to  cryogenic  systems  of  semiclassical  behavior  (Ba(Ar)„  clusters). 
Applications  to  Mg  in  solid  hydrogen,  both  partly  completed,  will  be  discussed  in  a  later 
section. 

2.1.2  The  CSP  Method  and  its  Extensions:  First  is  described  the  method  for  adiabatic 
systems.  The  method  was  outlined  in  (Refs.  2-4),  and  is  based  on  the  Classical  Separable 
Potential  (CSP)  approximation.  As  a  first  step  preceeding  the  quantum  calculation,  a  full 
classical  MD  simulation  of  the  same  process  was  carried  out.  Let  V(q{, . .  .,qN)  be  the 
full  potential  function  of  the  system,  and  let  q\a)(t), . . . ,  q\a\t)  be  the  trajectory  corre¬ 
sponding  to  the  a  initial  conditions.  For  each  mode  qh  an  effective  potential  for  that 
mode  is  defined: 

Vi{q»t)=t  . q^it),q^{t),...,q^\ t))w™ -V(t)  (1) 

a=  1 

The  coordinate-independent  quantity  V(r)  is  chosen  for  convenience  so  that  the  mean 
value  over  the  trajectories  of  the  sum  of  the  single-mode  potentials  at  time  t  equals  the 
value  of  the  full  potential  of  the  system  at  that  time.  w(ff)  in  Equation  (1)  is  the  weight  of 
the  a  initial  condition,  and  n  is  the  number  of  trajectories  used.  The  single-mode  effective 
potentials  are  used  to  generate  single-mode  wavefunctions,  or  time-dependent  "nuclear 
orbitals",  by  solving  the  time-dependent  Schrodinger  equation: 

ih  =  [Tj  +  Vj(qp  tWj(qj,  t )  (2) 

where  Tj  is  the  kinetic  energy  of  mode  j.  A  first,  but  often  very  useful  approximation  for 
the  total  wavefunction  in  terms  of  the  "nuclear  orbitals"  is  given  by  the  product: 

y'(qi,--',qN>t)  =  m>j(qj,t )  (3) 

j 

This  approximation  is  similar  in  some  respects  to  the  TDSCF  approximation  (Ref.  8)  and 
has  about  the  same  accuracy  -  but  is  incomparably  more  efficient  for  large  systems.  It  is 
readily  applicable  to  large  polyatomics  and  to  condensed  phases.  The  great  computational 
simplicity  of  the  method,  and  of  obtaining  the  "nuclear  orbitals",  together  with  other 
ideas,  makes  it  feasible  to  proceed  to  much  higher-level  approximations,  using  a  sum  of 
terms  of  the  product  type  Equation  (3): 

V(q , , . . . ,  qN ,  t)  =  X  Ck(t)  Tl/fCqj,  t )  (4) 

k  J=l 
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With  a  sufficient  number  of  terms,  this  "Configuration  Interaction"  expansion  in  terms  of 
the  time-dependent  nuclear  orbitals  can  be  exact  in  principle.  As  practically  applicable,  it 
can  provide  an  excellent  approximation  even  for  large  systems.  The  equations  of  motion 
for  the  Ck(t)  are  obtained  by  substitution  of  Equation  (4)  into  the  time-dependent 
Schrodinger  equation. 

The  basis  for  the  great  efficiency  of  the  CSP  is  in  the  very  efficient  production  of 
the  time-dependent  orbitals  for  atomic  motions,  which  in  turn  is  due  to  the  use  of  classi¬ 
cal  MD  to  construct  the  single-mode  effective  potentials,  from  which  the  orbitals  are 
obtained.  This  is  vastly  superior  to  the  effort  required  in  TDSCF,  where  effective  poten¬ 
tials  are  obtained  in  a  quantum  mechanical  self-consistent  procedure,  involving  computa¬ 
tionally  costly  integrations.  The  classical  MD  calculation  can  be  effective  also  in  suggest¬ 
ing  which  of  the  enormous  number  of  possible  configurations  in  Equation  (4)  should  be 
retained.  We  do  not  give  the  details  here,  but  a  quantity  that  plays  a  key  role  is  the  "corre¬ 
lation  potential" 

AV(qu  •  •  • ,  qN,  t)  =  V(qx ,  ■  ■  ■ ,  qN)  -  X  ?«(*,•.  0  (5) 

i 

where  V  is  the  full  potential  of  the  system,  and  the  V^t)  are  the  CSP  single-mode  effec¬ 
tive  potentials.  AV  represents  the  dynamical  correlation  between  the  various  degrees  of 
freedom,  not  present  in  the  simple  CSP.  A  set  of  orbitals  of  Equation  (4)  that  provides  a 
sufficient  basis  to  represent  AV  is  enough  for  an  accurate  solution  of  the  Schrodinger 
equation.  Calculations  in  the  "Configuration  Interaction"  extension  of  CSP  can  practically 
include  several  hundreds  of  terms  in  Equation  (4),  and  this  is  by  the  evidence  of  numeri¬ 
cal  tests,  sufficient  for  quite  good  accuracy  even  in  large  systems. 

For  some  processes  and  systems,  especially  when  the  required  time  scales  are  very 
short  (say,  t  «  1000  femtoseconds),  the  simplest  level  of  CSP  suffices.  However,  for 
processes  involving  nonadiabatic  transitions  between  different  electronic  states,  this  is 
never  the  case.  The  method  then  requires  important  generalizations: 

(i)  Rather  than  use  classical  MD,  one  determines  trajectories  from  a  semiclassical  MD 
that  allows  for  transitions  between  different  electronic  states.  We  used  for  this  purpose 
Tully’s  semiclassical  surface-hopping  method.(Ref.  9)  From  these  trajectories  one  can 
construct  effective  single-mode  Hamiltonians,  and  use  the  latter  to  construct  the  time- 
dependent  single-mode  orbitals. 

(ii)  The  effective  Hamiltonians  constructed  in  this  case  already  include  the  effects  of 
nonadiabatic  transitions.  The  full  approach  requires  in  this  case  the  use  of  a  wavefunction 
that  has  at  least  n  terms,  each  of  which  is  a  simple  product  of  orbitals  (single-mode 
wavepackets),  n  being  the  number  of  electronic  states  involved.  The  form  of  the 
wavefunction  in  the  nonadiabatic  case  is,  in  the  simplest  version: 

¥(r, ql,---,qN,t)='£  Ck{t)y/k( r, q) %(lc\q, t)  (6) 

k=i 


where 

Zm(q,t)=n/f’(qj,0  (7) 

In  Equation  (6),  r  denotes  the  electronic  coordinates  and  q  the  nuclear  ones.  The  Ck(t) 
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are  time  dependent  coefficients,  and  are  determined  by  solution  of  the  Schrodinger  equa¬ 
tion.  The  orbitals  <pf(qj,t)  are  determined  by  solving  the  time-dependent  equation  for 
the  effective  Hamiltonians  associated  with  the  individual  modes  qj.  The  foil  procedure  in 
the  nonadiabatic  case  is  described  in  Ref.  4.  Tests  for  small  model  systems  have  shown 
that  this  is  an  approach  of  very  good  accuracy. 

2.1.3  Applications:  The  CSP  method  and  its  extensions  were  tested  extensively  for 
small  model  systems,  and  were  applied  already  to  several  realistic  cases,  all  in  the  frame¬ 
work  of  the  present  project.  The  earliest  tests  of  the  simple  SCP  were  for  photoinduced 
processes  in  collinear  model  clusters,  e.g.  photodetachment  dynamics  of  I_(Ar)2,  and 
photodissociation  of  HX  in  Rg  •  •  •  HX.(Ref.  2)  Many  of  the  more  recent  tests  are  also  for 
triatomic  systems,  since  for  such  small  systems  there  are  exact  reference  calculations  by 
other  time-dependent  methods.  Estimates  of  the  errors  involved  can  be  carried  out  in  the 
framework  of  the  method  itself  (based  on  calculations  of  correction  terms),  and  these 
were  used  for  some  of  the  realistic  applications. 

The  first  applications  to  realistic  systems  were  for  processes  in  myogenic  clusters 
of  Ne  and  Ar  atoms,  rather  than  for  hydrogen  systems.  The  less  pronounced  nature  of  the 
quantum  effects  in  the  "semiclassical"  clusters  made  applications  and  interpretations  with 
the  new  method  simpler,  since  the  physics  of  the  processes  involved  is  easier  to  under¬ 
stand.  Discussed  here  is  such  an  application  to  a  heavy-atom  cryogenic  cluster,  and  defer 
to  Section  (2.8)  results  on  simulations  for  Mg  in  H2  clusters,  which  is  a  system  of  direct 
HEDM  relevance.  The  application  described  here  is  of  the  photoexcitation  dynamics  of 
Ba(Ar)„  clusters,  with  10  <  n  <  20.  A  major  reason  for  choosing  this  system  is  that  one 
expects  it  to  be,  at  least  to  some  extent,  a  useful  (if  more  classical)  model  for  Mg(H2)„. 
Another  important  reason  is  that  interaction  potentials  of  good  quality  were  available  in 
this  case  from  empirical  sources.(Ref.  10)  Finally,  for  large  Ba(Ar)„  clusters,  experimen¬ 
tal  results  were  available  from  the  group  of  Mestdagh,  Visticot  et  al.{ Refs.  10,11)  in 
Saclay,  France.  In  fact,  in  a  recent  joint  study  with  that  group,  a  nearly  quantitative  inter¬ 
pretation  for  their  findings  was  provided.(lO) 

The  application  is  described  in  detail  in  a  paper  by  Jungwirth  and  Gerber.(Ref.  4) 
Using  the  new  method  described  above,  the  Ba(S)  Ba(P)  excitation  in  Ba(Ar)„  clus¬ 
ters,  for  n  —  10,20  was  simulated.  The  structure  of  the  clusters  in  the  electronic  ground 
state  roughly  corresponds  to  a  Ba  atom  adsorbed  on  a  surface  of  Ar  atoms.  For  this  struc¬ 
ture  the  3  excited  Retype  singlet  state  can  be  classified  as  one  I-type  state  (where  the  ori¬ 
entation  of  the  p-orbital  of  the  Ba  is  roughly  normal  to  the  surface)  and  two  n-type  states 
(where  the  ^-orbital  density  is  approximately  parallel  to  the  surface).  In  the  dynamical 
calculations  excitation  into  the  highest  (E-type)  state  was  considered,  and  followed  the 
development  of  the  system  in  time.  The  CSP  results  were  compared  with  those  of  semi¬ 
classical  surface  hopping  calculations  (Tully’s  method).  The  nonadiabatic  version  of  the 
CSP,  with  3  configurations  in  the  wavefunction  expansion  (Equation  (6))  was  used.  Tests 
indicate  that  the  CSP  results  are  very  accurate  in  this  case,  and  the  validity  of  the  method 
holds  for  times  up  to  f— 1. 5  fs  and  probably  longer.  Some  of  the  main  results  are  shown  in 
Figs.  (1),  (2).  To  emphasize  some  of  the  key  points: 

(i)  The  very  symmetric,  well-defined  structure  of  the  initial  state  is  destroyed  by  the 
excited  state  motions  within  a  timescale  of  t  <  500  fs.  As  Figure  1  for  the  Ba(Ar)10 
shows,  the  initial  state  corresponds  to  3  very  narrow  peaks  in  the  Ba-Ar  distance 


7 


Ba-Ar  ^Ba-Ar 


t  =  Qfs 


r  (Angstrom) 


r  (Angstrom) 


t  =  1 000  fs 

2 


1 

0 

0 

5  6  7  8  9  10  11  12  13 
r  (Angstrom) 


5  6  7  8  9  10  11  12  13 
r  (Angstrom) 


Figure  1 

Time  Evolution  of  the  Ba-Ar  Atom- Atom  Distribution  Function  for  the  Ba(Ar)10  Cluster 
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distribution,  reflecting  a  well-defined  relatively  rigid  cluster  structure  where  the  Ba  is 
adsorbed  on  top  of  an  approximately  planar  (Ar)10.  The  3  main  Ba-Ar  distances  are  only 
faintly  recognizable  at  t  =  500  fs,  and  the  whole  structural  distribution  has  become 
smeared.  The  trend  continues  for  longer  times  excitation,  and  for  t  =  1500  fs  the  distribu¬ 
tion  is  so  broad  and  smeared  that  it  is  meaningless  to  attribute  a  structure  to  the  cluster. 
Most  of  the  smearing  is  due  to  the  motion  of  the  Ba  atom  away  from  the  (Ar)10  surface 
(the  probability  for  direct  dissociation  of  the  cluster  in  this  process  is  about  80%),  but 
there  are  also  changes  in  the  structure  of  the  (Ar)10  "surface". 

(2)  Figure  2  shows  the  relative  populations  (branching  ratios)  of  the  three  different  elec¬ 
tronic  P-states  of  the  Ba(Ar)10  system,  as  a  function  of  time.  The  figure  shows  that  elec¬ 
tronic  relaxation  from  the  £  to  the  IT  states  becomes  significant  after  a  timescale  of  -  1 
picosecond  or  less,  the  released  energy  going,  of  course,  into  the  vibrational  degrees  of 
freedom  of  the  cluster.  The  relaxation  is  dominantly  into  the  state  intermediate  in  energy  - 
the  population  of  the  lowest  O  state  remains  extremely  small  at  least  through  the  first 
1500  fs  of  the  process. 

(3)  For  Ba(Ar)10,  direct  dissociation  is  the  dominant  process  (-80%).  Only  about  20%  of 
the  wavepackets  describe  electronic  relaxation  into  the  11  state  and  stabilization  of  the 
cluster.  In  the  case  of  Ba(Ar)20,  electronic  relaxation  greatly  dominates  over  direct  disso¬ 
ciation.  (Direct  dissociation  is  about  20%  only).  This  trend  becomes  stronger  for  larger 
Ba(Ar)„  clusters:  For  the  experimental  clusters  of  Visticot  et  al.  (Refs.  10,11),  with 
n  »  100,  our  simulations  showed  that  evaporation  of  Ba  from  the  cluster  after  excitation 
is  a  fairly  rare  event,  with  a  probability  of  several  percent  at  most.  (Ref.  12)  As  cluster 
size  increases,  so  does  the  binding  energy,  and  therefore  the  decrease  of  direct  dissocia¬ 
tion.  Electronic  relaxation  becomes  slower  also  with  cluster  size,  but  less  so  than  direct 
dissociation.  It  has  a  timescale  of  103  fs  only  for  the  very  large  clusters.  These  conclu¬ 
sions  on  electronic  relaxation  of  excited  metals  in  cryogenic  media  are  of  importance  for 
HEDM  purposes,  since  electronic  excitation  has  been  used  to  characterize  and  study 
metal-doped  solid  hydrogen  systems.  (Refs.  13, 14) 

(4)  Figure  2  shows  that  the  semiclassical  surface-hopping  method  can  be  in  significant 
error  in  predicting  electronic  branching  ratios.  Qualitatively,  the  prediction  of  this  method 
of  the  populations  in  time  looks  reasonable,  but  the  result  for  the  population  of  the  inter¬ 
mediate  level  is  off  by  a  factor  of  2  for  t  =  1. 5  fs.  (The  error  for  the  lowest  state  is  even 
greater,  due  to  the  fact  that  one  deals  with  a  very  low  population  there,  the  magnitude  of 
which  is  hard  to  get  accurately).  On  the  other  hand,  the  semiclassical  method  is  success¬ 
ful  in  predicting  quite  reasonably  the  change  of  the  structure  of  the  cluster  in  time  -  see 
Figure  1.  Thus,  even  for  heavy  atom  clusters,  the  semiclassical  surface  hopping  method 
must  be  used  with  caution,  and  only  for  certain  properties.  Whenever  possible  in  doubtful 
cases,  use  of  the  new  method  proposed  here  is  recommended. 

2.1.4  State  of  the  Method  and  Future  Directions:  The  algorithm  and  code  for  the  sim¬ 
plest  level  CSP  have  already  been  developed  to  a  level  of  an  efficient,  user-friendly  tool. 
Simulations  for  several  systems  with  over  150  coupled  degrees  of  freedom  have  been  car¬ 
ried  out.  Code  development,  with  emphasis  on  implementation  on  parallel  computers, 
continues.  Anticipated  are  applications  for  systems  with  over  500  atoms  (1500  degrees 
of  freedom)  in  a  matter  of  months  at  most.  When  simple  CSP  suffices,  it  is  already  a  very 
powerful  practical  tool  for  Quantum  Molecular  Dynamics.  The  extension  of  the  method 
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to  include  many  configurations,  as  in  Equation  (4),  has  the  advantage  of  making  the 
approach  very  accurate  (in  fact,  with  a  sufficient  number  of  terms  in  Equation  (4),  it 
should  be  "numerically  exact").  The  practical  performance  in  accuracy  of  the  method,  in 
all  test  cases,  is  very  promising  and  exciting.  However,  the  largest  calculations  with  it  are 
still  confined  to  ~  15  atoms  (less  than  50  degrees  of  freedom).  Again,  developments  on 
this  are  very  fast,  and  systems  of  the  order  of  a  few  hundreds  of  modes  should  be  devel¬ 
oped  within  a  year  or  so.  Presently,  there  is  no  other  method  of  similar  accuracy  for  the 
system  sizes  treated. 

2.2  Quantum  and  Semiclassical  TDSCF  Simulations  of  Cryogenic  Systems 

The  original  plan  for  this  research  was  to  develop  a  time-dependent  quantum 
method  for  simulations  of  many-particle  systems  based  on  the  Time-Dependent  Self- 
Consistent  Field  (TDSCF)  approximation.  (Ref.  8)  A  CSP  was  developed  as  the  simula¬ 
tion  tool  for  the  purposes  of  this  project,  and  as  a  powerful  and  practical  algorithm  for 
time-dependent  quantum  simulations  of  large  systems  in  general.  Nevertheless,  much  of 
the  research  on  this  project  dealt  with  TDSCF,  and  some  very  useful  results  were 
obtained  with  this  method.  As  seen  below,  the  method  is  useful  in  some  practical  applica¬ 
tions,  including  certain  specific  systems  of  HEDM  interest,  although  in  general  it  does 
not  seem  to  compete  favorably  with  CSP.  Presented  now  are  the  TDSCF  results  obtained 
in  the  project,  with  some  comments  and  conclusions  on  HEDM  applicability.  The  TDSCF 
approximation  has  been  known  in  general  formulation  from  the  early  stages  of  quantum 
theory,  and  has  been  in  extensive  use  by  several  groups,  including  the  author’s  group. 
This  refers  to  a  review  by  Ratner  and  Gerber  on  the  principles  of  the  method.  (Ref.  8) 
Examined  here  are  conclusions  and  findings  obtained  in  this  project. 

As  seen  below,  simulations  have  shown  that  TDSCF  is  of  good  accuracy  for  very 
short  timescale  dynamics  (t  «  lOOOfs)  even  for  hydrogen  clusters.  (Ref.  15)  For  longer 
times,  extensions  are  required  if  accuracy  of  the  description  is  to  be  maintained.  The 
TDSCF,  as  the  simplest  version  of  the  CSP,  represents  the  total  wavefunction  of  the  sys¬ 
tem  as  a  product  of  time-dependent  single-mode  wavepackets  (nuclear  orbitals).  Exten¬ 
sions  of  TDSCF  have  been  proposed  that  require  sums  of  terms  in  the  wavefunction 
expression  with  each  term  being  a  product  of  single-mode  wavepackets  (Refs.  16-20), 
such  as  multiconfiguration  TDSCF.  (Refs.  18-20)  The  failure  of  TDSCF  is  only  a  failure 
of  computational  efficiency: 

(1)  For  large  systems  and  for  complicated  interaction  potentials,  the  TDSCF  approxima¬ 
tion  itself  becomes  computationally  very  costly.  The  calculation  of  single-mode  effective 
potentials  involves,  for  realistic  interactions  between  the  particles,  the  evaluation  at  each 
time  point  of  multidimensional  integrals.  As  a  result,  even  for  systems  with  a  modest 
number  of  modes,  say  -10,  the  CSP  is  by  an  order  of  magnitude  computationally  superior 
to  TDSCF. 

(2)  In  the  multiconfiguration  version  of  the  theories,  the  handicap  of  TDSCF  is  even  more 
severe.  Multiconfiguration  TDSCF  is  a  very  accurate  high  level  algorithm  for  molecular 
dynamics,  and  works  very  well  also  for  HEDM-type  systems,  but  is  essentially  limited  to 
systems  of  several  atoms  (except  if  one  has  the  advantage  of  very  special,  simple  poten¬ 
tials  which  are  usually  unrealistic)  CSP  has  an  enormous  advantage  in  the  efficiency  of 
generating  the  single-mode  wavepackets,  and  also  in  assessing  how  many  and  which 
types  of  configurations  should  be  included  (these  are  due  to  the  use  of  classical  MD  in 
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generating  effective  one-mode  potentials,  hence  orbitals,  and  in  providing  guidance  for 
the  quantum  calculation.)  TDSCF  becomes  thus  a  practical  method  either  for  simple  sys¬ 
tems,  or  when  special  simplifying  potentials  can  be  used. 

Validity  of  the  TDSCF  approximation  for  the  dynamics  of  hydrogen  clusters:  In  order  to 
test  whether  TDSCF  holds  for  the  dynamics  of  HEDM-type  materials,  we  carried  out 
comparisons  between  TDSCF  and  exact  calculations  for  small  metal  atom  -  hydrogen 
clusters.  (Ref.  15)  One  of  the  systems  for  which  the  calculations  were  done  is  a  collinear 
model  of  the  cluster  B(H2)2,  having  the  structure  H2  -  B  -  H2.  For  such  a  small  system, 
numerically  exact  solutions  of  the  time-dependent  Schrodinger  equation  is  possible.  The 
TDSCF  approximation,  the  accuracy  of  which  depends  on  the  choice  of  coordinates  since 
the  approximations  involve  a  mutual  separation  of  the  modes,  was  applied  in  normal 
coordinates.  In  the  calculations,  a  non-stationary,  but  bound,  initial  state  of  the  cluster 
1  >  =  0)  was  propagated  in  time  exactly,  and  also  in  the  TDSCF  approximation. 

The  autocorrelation  function  of  the  wavefunction  <  'i'(qi,q2,t)\'¥(qi,q2,0)  >  was  one  of 
the  quantities  used  to  compare  the  TDSCF  and  the  exact  results.  This  is  a  phase-sensitive 
quantity,  and  quite  demanding  for  the  approximation.  It  was  found  that  the  TDSCF 
approximation  remains  valid  within  reasonable  accuracy  to  times  of  the  order  of  t  <  50fs 
for  B(H2)2.  In  this  domain,  the  autocorrelation  function  computed  from  TDSCF  is  accu¬ 
rate  to  80%  or  better.  For  longer  times,  the  TDSCF  approximation  error  builds  up,  and 
application  will  therefore  require  extensions  of  TDSCF,  as  discussed  above.  These  calcu¬ 
lations  of  Li  and  Gerber  thus  show  that  the  TDSCF  is  applicable  for  HEDM-type  sys¬ 
tems,  in  ultrafast  time  domains.  As  mentioned  earlier,  the  problem  in  applying  the 
TDSCF  is  the  computational  effort  that  is  involved  for  large,  real  systems.  TDSCF  calcu¬ 
lations  are  feasible  and  efficient  for  large  collinear  model  clusters,  and  such  calculations 
were  carried  out  by  Li  and  Gerber  for  (Ne)„,  n  <  16.  (Ref.  21)  For  large  clusters  in  3D 
geometry,  the  calculations  became  very  costly,  and  TDSCF  is  computationally  much  infe¬ 
rior  to  CSP,  except  when  certain  specific,  advantageous  potential  functions  can  be  used, 
as  discussed  below. 

TDSCF  simulations  for  large  systems  with  quartic  force  fields:  In  cases  when  the  multidi¬ 
mensional  potential  surface  can  be  expanded  as  a  power  series  in  the  coordinates  with  a 
sufficiently  low  maximal  power,  the  TDSCF  was  found  to  be  practically  applicable  and 
efficient,  and  can  be  used  for  Quantum  Molecular  Dynamics  simulations  of  large  systems 
in  3D  geometry.  It  was  discovered  that  for  the  clusters  (Ar),3  and  (Ne)I3  for  instance,  a 
quartic  representation  of  the  potential  in  the  normal  modes  is  of  sufficient  accuracy  with 
low  energies.(20)  Thus,  a  potential  function  of  the  form: 

V(qu-  ■■,qN)=  X  Cn(l),  •  •  • ,  n(N)(qi)n(l)  •  •  •  (qN)n{N)  (8) 

n(i)>0 

was  used  where  the  C„(1),  •  •  •  ,n(N)  are  the  expansion  coefficients,  and 
n(l)  +  ---  +  n(iV)  <4.  Using  this  force  field,  simulations  of  (Ar)13,(Afe)13  were  carried 
out  and  other  large  clusters  in  full  3D  geometry.  (Ref.  22)  There  were  semiclassical 
TDSCF  simulations,  where  in  addition  to  the  TDSCF  approximation,  a  semiclassical 
Gaussian  approximation  is  used  for  the  single  mode  wavepackets.  The  semiclassical 
approximation  holds  in  cases  when  the  anharmonicity  is  not  very  large.  To  test  the 
method  for  a  large,  realistic  system  equilibrium  simulations  were  carried  out,  that  is,  the 
semiclassical  TDSCF  for,  e.g.,  (Ar)13,(Ne)13  was  run  until  the  system  reached  equilib¬ 
rium,  for  which  criteria  was  provided.  (Ref.  22)  The  results  were  compared  with  the 
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Feynman  Path  Integral  (FPI)  simulations  which  was  carried  out  for  the  same  system:  The 
FPI  is  "numerically  exact",  but  cannot  be  used  in  practice  for  real-time  simulations,  only 
for  equilibrium  ones  (for  realistic  cases).  The  FPI  could  be  used  to  test  the  semiclassical 
TDSCF.  Figure  3  shows  a  comparison  between  the  radical  distribution  function  g(r) 
computed  from  TDSCF,  and  the  corresponding  result  from  FPI,  (20)  for  (Ne)13.  Both  cal¬ 
culations  are  for  T  =  2K.  The  agreement  is  excellent,  establishing  the  validity  of  semi¬ 
classical  TDSCF  for  this  system.  For  hydrogen  clusters  and  solids  under  normal  condi¬ 
tions,  the  anharmonicity  is  so  high  that  expansion  (Equation  (8))  takes  many  terms  to 
converge,  and  the  approach  is  impractical.  However,  for  hydrogen  at  high  pressure,  P  > 
5kbar,  say,  the  anharmonicity  is  no  longer  so  severe,  and  expansion  (8)  can  be  used  with  a 
reasonable  number  of  terms.  Thus,  it  is  shown  that  TDSCF  can  be  applied  for  time- 
dependent  quantum  simulations  of  solid  hydrogen  and  doped  hydrogen  materials  at  high 
pressure.  The  method  can  also  be  applied  to  other  cryogenic  systems,  like  solid  Ne,  that 
are  not  of  HEDM  applications  but  are  useful  as  models  in  HEDM  research.  In  conclusion, 
it  is  shown  that  the  TDSCF  method  can  be  useful  as  a  Quantum  Molecular  Dynamics 
simulation  tool  in  HEDM  research,  but  only  for  a  restricted  class  of  systems. 

2.3  Algorithms  for  Equilibrium  Simulations  of  Quantum  Clusters 

The  tasks  set  at  the  initiation  of  this  research  project  also  required  simulations  of 
ground-states  or  thermal  equilibrium  states  of  HEDM-related  systems.  For  instance,  cal¬ 
culations  of  time-dependent  dynamical  processes  in  the  systems  of  interest  often  require 
simulations  of  equilibrium  states  as  the  appropriate  initial  states  for  the  dynamical  pro¬ 
cesses.  Essentially,  the  photoinduced  processes  considered  in  this  work  begin  from  an 
equilibrium  system.  Simulations  for  two  types  of  equilibrium  states  were  carried  out:  (1) 
The  vibrational  ground  state  of  HEDM-related  clusters  or  solids;  (2)  Thermal  equilibrium 
states  at  finite  temperature  of  these  systems.  For  both  types  of  calculations  methods  well- 
established  in  the  literature  were  used,  but  developed  in  each  case  its  own  specific  algo¬ 
rithms  and  codes,  which  have  specific  advantages  for  our  purposes  compared  with  other, 
exciting  ones  based  on  the  same  basic  methods.  The  simulations  using  these  codes  played 
a  major  role  in  accomplishing  the  tasks  as  set  in  projects  detailed  in  the  Program  Plan.  In 
addition,  the  availability  of  the  codes  may  prove  a  useful  future  resource  for  HEDM 
research  by  other  groups. 

2.3.1  Diffusion  Quantum  Monte  Carlo  (DQMC)  Simulations  of  Vibrational 
Ground  States  of  HEDM  Systems:  DQMC,  and  related  methods  such  as  Variational 
Quantum  Monte  Carlo,  are  amongst  the  powerful  available  techniques  for  calculations  of 
the  ground  state  of  large  quantum  systems.  (Refs.  23-28)  Our  work  (in  cooperation  with 
M.  Alexander)  was  the  first  to  extend  the  method  for  systems  containing  an  open-shell,  p- 
state  atoms.  (Ref.  29)  Also,  a  code  was  developed  for  implementing  the  work  on  a  mas¬ 
sively-parallelized  computer  (Ref.  30)  and  indeed  the  work  was  mostly  earned  out  on  the 
MasPar2  at  U.C.  Irvine,  a  system  of  4096  nodes.  This  resulted  in  a  capability  to  simulate 
large  systems  very  efficiently.  In  the  simulations  described  in  the  forthcoming  sections 
DQMC  was  used  to  obtain  structural  and  energetic  properties  of  potential  HEDM  sys¬ 
tems.  Table  1  presents  data  on  the  energies  of  HEDM-related  clusters  obtained  in  our 
calculations,  and  more  extensive  data  and  analysis  of  the  main  systems  is  presented  later 
on  in  this  report. 
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Weigth  of  distance  R  [arbitrary  unit] 


Neon  -  Neon  Distance  R  in  [Angstrom] 


Figure  3 

The  Ne-Ne  Distance  Distribution  function  for  (Ne)i3  at  T=1K.  Semiclassical  Molecular 
Dynamics  Results  are  Compared  with  Feynman  Path  Integral  Simulations. 

From  Ref.  (22) 
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Table  1.  DQMC  Results  for  the  Cohesion  Energy  of  Quantum  Clusters 


Cluster 

Vibrational  Energy 
A)(cm_1) 

h2-h2 

1.86 

(H2)3 

6.867 

(H2)4 

14.06 

(H2)5 

25.01 

(h2)6 

34.65 

(H2)io 

86.99 

(H2),3 

134.5 

(H2)45 

515.54 

(H2)ioo 

843.02 

Hg-H2 

13.65 

Hg  ~  (Hq)2 

31.80 

Hg  — (H2)io 

212.3 

Hg  -  (H2)12 

282.8 

Hg  —  (H2)24 

485.0 

Hg  -(H^o 

687.7 

(D  2)3 

16.0 

Li-H2 

8.5 

Li  -  (H2)2 

3.6 

o-h2 

15.00 

o-d2 

19.80 

f-h2 

14.06 

(At)  i3 

4586.8 
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2.3.2  Feynman  Path  Integral  Simulations  of  Quantum  Systems  at  Thermal  Equi¬ 
librium:  Extensive  work  on  this  topic  was  done  for  Li(H2)„  clusters  by  Klein  and 
coworkers.  (Ref.  31)  The  FPI  simulations  were  applied  to  structural  and  energetic  proper¬ 
ties  of  Mg(H2)„  and  other  quantum  clusters.  (Ref.  32)  The  algorithm  (Ref.  22)  followed 
the  approach  of  References  (33)  and  (34).  Also,  the  algorithm  was  recently  extended  for 
FPI  simulations  of  systems  containing  open-shell  atoms.  Furthermore,  the  algorithm  was 
implemented  in  codes  for  distributed  computers.  The  FPI  simulation  data  are  given  in 
subsequent  sections. 

2.4  Structural  Properties,  Stability  and  Dynamics  of  B(H2)„  Clusters 

The  results  obtained  in  this  part  are  amongst  the  most  important  ones  of  the 
research  project,  since  they  support  the  suggestion  that  boron  in  hydrogen  may  be  one  of 
most  promising  candidates  for  a  HEDM  propellant.  Indeed,  the  results  seem  most  encour¬ 
aging.  This  work  was  done  in  cooperation  with  M.  Alexander  and  his  group.  The  B  -  H2 
interaction  potential  used  was  taken  from  the  ab  initio  calculations  of  Alexander.  (Ref. 
35)  The  calculations  were  done  by  Diffusion  Quantum  Monte  Carlo  (DQMC),  adapted  to 
deal  with  a  system  that  includes  an  open-shell  P-state  atom.  An  extensive  account  of  the 
results  is  given  in  Reference  (29),  and  only  a  brief  summary  of  the  findings  is  given  here. 
The  calculations  of  Reference  (29)  were  for  B(H2)„,  with  1  <  n  <  8.  Calculations  for  clus¬ 
ters  with  a  much  larger  number  of  H2  molecules  (up  to  n  =  100)  are  under  way.  The  pre¬ 
liminary  results  support  the  conclusions  below,  and  are  based  on  the  small-cluster  data. 
Listed  below  is  a  summary  of  the  findings,  and  of  the  conclusions  regarding  structural 
properties,  energetic  stability,  boron  atom  diffusitivity  and  vibrational  properties  of 
boron-hydrogen  clusters: 

(1)  B(H2)„  clusters  are  energetically  very  stable,  compared  with  pure  H2  clusters  or  with 
clusters  of  5-state  atoms  (e.g.,  Li,  Mg)  with  hydrogen  molecules.  As  Table  2  shows,  one 
expects  B(H2)„  (at  least  for  n<6),  to  remain  stable  well  above  the  melting  point  of  pure 
solid  hydrogen.  The  low  energy  of  the  minima  of  the  potential  functions  are  offset  in  part 
by  the  high  zero-point  energies.  Nevertheless,  on  the  basis  of  the  results,  B  atoms  in  solid 
H2  should  undergo  strong  complexation  with  neighboring  hydrogen  molecules.  At  a  tem¬ 
perature  of  4K,  considered  in  some  HEDM  designs,  the  solvated  boron  atom  should  be 
extremely  stable. 

(2)  To  further  quantify  the  energetic  stability  of  B(H2)n  clusters,  and  the  effective  solva¬ 
tion  of  boron  by  hydrogen,  the  accretion  energy  and  the  solvation  energy  of  the  B(H2)n 
clusters  was  considered.  The  accretion  energy  is  defined  as  the  energy  gained  by  adding 
one  more  H2  molecule  to  B(H2)„.  The  solvation  energy  is  defined  as  the  energy  gained  by 
adding  a  B  atom  to  a  pure  H2  cluster.  The  accretion  and  solvation  energies,  as  a  function 
of  cluster  size,  are  shown  in  Figure  4.  The  special  stability  of  the  n  =  5  cluster  is  evident 
from  the  results.  The  values  of  accretion  and  of  solvation  energy  for  all  the  clusters 
shown  are  very  large  on  the  scales  of  such  quantities  found  for  M(H2)„  clusters,  where  M 
is  an  alkali  or  an  earth-alkali  metal  (Li  and  Mg  were  considered). 
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Table  2.  Calculated  Potential  Energy  Minima 
and  Lowest  Vibrational  Energies  of  B(H2)„ 


n 

Potential  Minimum 

Vibrational  Energy 

i 

-61.3 

-25.4 

2 

-146.8 

-47.5+1.0 

3 

-233.7 

-74. 8±2. 5 

4 

-309.0 

-97.4±1. 9 

5 

-410.8 

-137.9+3.5 

6 

-523.7 

-156. 4±4. 3 

(Energies  are  given  in  cm  ’) 

(3)  The  B(H2)„  clusters  are  very  anisotropic  in  structure.  This  stems  from  the  dependence 
of  the  boron-hydrogen  interaction  on  the  unpaired  p-orbital  of  the  B  atom:  The  interac¬ 
tion  is  considerably  deeper  when  the  H2-B  axis  is  normal  to  the  direction  of  the  p-orbital 
(7r-state  in  a  pseudo-atom  treatment  of  H2).  Therefore,  for  low  n,  the  lowest  energy  struc¬ 
tures  of  B(H2)„  are  planar,  with  the  plane  of  the  H2  molecules  normal  to  the  p-orbital.  It 
must  be  kept  in  mind,  however,  that  the  structure  of  none  of  the  clusters  is  rigid,  but 
rather  fairly  delocalized.  The  structures  corresponding  to  minima  of  the  potential  surfaces 
are  shown  in  Figure  5  for  n  =  5, 6, 8.  In  each  case  are  shown  the  structure  corresponding 
to  the  global  minimum,  and  the  structure  corresponding  to  local  low-lying  minima.  The 
B-atom  is  at  the  center  of  the  cluster  for  the  global  minimum  of  B(H2)6,  but  not  in  the 
other  cases.  The  structures  corresponding  to  the  local  minima  of  the  potential  surface  dif¬ 
fer  in  some  cases  substantially  from  the  global  minimum  structure.  Not  in  all  cases  is  this 
important  since  the  vibrational  wavefunctions  are  very  delocalized,  and  can  be  spread 
over  several  minima.  At  least  some  of  the  local  minima  are,  however,  outside  the  delocal¬ 
ization  range,  and  thus  the  existence  of  multiple  structural  "isomers"  is  an  important 
property  of  B(H2)„.  In  conclusion,  one  important  finding  are  the  highly  nonspherical, 
anisotropic  structures  of  B(H2)„  and  another  is  the  fact  that  above  the  lowest  temperatures 
several  different  structures  may  be  populated. 

(4)  Data  on  the  zero-point  quantum  fluctuation  of  the  structures  of  B(H2)„,  is  presented  in 
Figure  6.  The  upper  part  shows  the  distribution  of  B  atom  distances  from  the  center  of  the 
cluster.  The  middle  part  of  Figure  6  shows  the  distribution  of  H2  distances  from  the  cen¬ 
ter.  Finally,  the  lower  part  shows  the  angular  distribution  of  H2  molecules  with  respect  to 
the  direction  of  the  unpaired  p-orbital  of  the  B.  It  is  evident  that:  (a)  The  structural  distri¬ 
butions,  due  to  the  zero-point  quantum  fluctuations  are  large  (although  the  structures  are 
less  delocalized  than  for  the  corresponding  pure  (H2)„  clusters);  (b)  The  B  atom  is  much 
closer  to  the  center  of  the  cluster  than  the  H2  atoms  that  "surround"  it;  (c)  The  zero-point 
motions  are  not  sufficient  to  "smear  out"  the  anisotropic  structure  of  the  B(H2)„  clusters. 
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The  Structures  of  Global  and  Local  Minima  of  Clusters  B(H2)n.  The  Notation  n.m  is  Used  to 
Designate  Minimum  m  of  the  Cluster  with  n  Hydrogens,  m  =  0  is  the  Global  Minimum,  m  =  1  is 

the  Lowest  Local  Minimum,  etc. 
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Figure  6 

Structural  Distribution  for  B(H2)n  from  DQMC 
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(5)  The  role  and  extent  of  H2  orientation  within  the  clusters  were  determined  in  the  case 
of  B(H2).  While  in  B(H2)  there  is  more  anisotropy  of  the  H2  orientation  than,  e.g.,  in 
(H2)„,  or  in  Li(H2)„,  the  anisotropy  is  slight,  and  one  can  treat  the  H2  as  a  spherical  parti¬ 
cle  to  good  approximation. 

(6)  From  the  strong  binding  of  the  small  clusters  B(H2)n),  it  is  inferred  that  diffusivity  of 
B  atoms  in  solid  H2  will  be  extremely  low,  up  to  the  melting  point  of  H2:  The  stable  local 
B(H2)n  complexes  in  B-doped  solid  hydrogen  will  prevent  diffusion  of  the  boron.  This 
implies  that  solid  hydrogen  doped  with  B  should  be  stable  against  diffusion,  and  remain 
non-segregated,  preferably  up  to  the  melting  temperature  of  hydrogen.  This  is  an  impor¬ 
tant,  favorable  property  of  B  in  solid  hydrogen  for  HEDM  purposes.  This  part  of  the  pro¬ 
ject  accomplished  Tasks  1  and  4  of  the  Program  Plan  of  the  research  proposal. 

2.5  Structural  Properties,  Stability  and  Vibrational  Dynamics  of  Mg(H2)„  and 
Related  Clusters 

The  calculations  on  Mg(H2)„  clusters  were  motivated  by  the  estimates  of  P.  Carrick 
of  the  Phillips  Laboratory,  suggesting  that  Mg  in  solid  hydrogen  has  a  high  Specific 
Impulse,  and  should  be  of  potential  interest  as  a  HEDM  propellant.  The  objective  of  the 
calculations  we  carried  out  was  to  determine  the  energetic  stability,  the  structural  proper¬ 
ties  and  the  vibrational  properties  of  Mg(H2)n.  Calculations  on  the  structural  and  ener¬ 
getic  properties  of  metal-hydrogen  systems  were  previously  carried  out  by  M.  Klein  and 
coworkers  for  Li(H2)„  and  also  for  Li  in  solid  hydrogen.  (Refs.  31,36)  These  calculations, 
using  the  Feynman  Path  Integral  method  for  thermal  equilibrium  simulations  showed  that 
in  Li(H2)„  the  Li  is  very  weakly  adsorbed  on  the  surface  of  an  essentially  spherical  (H2)„. 
An  elegant  experimental  approach  has  made  it  possible  for  Fajardo  and  coworkers  to 
"implant"  Li  impurities  in  solid  hydrogen  and  study  their  spectroscopy.  (Refs.  13,14) 
Nevertheless,  the  results  of  Klein  and  coworkers  suggest  very  unfavorable  prospects  for 
preparing  appreciable  concentrations  of  isolated  Li  atoms  in  in  solid  hydrogen.  Prior  to 
obtaining  an  interaction  potential  for  Mg  -  H2,  (Ref.  37)  simulations  were  carried  out  for 
Hg(H2)„,  and  is  expected  to  prove  similar  in  properties  to  Mg(H2)„.  Although  Hg(H2)„  is 
not  a  HEDM  candidate,  its  properties  are  useful  in  throwing  light  on  the  expected  behav¬ 
ior  of  Mg(H2)„. 

(1)  Mg-H2  Interaction:  In  our  DQMC  simulations  on  Mg(H2)„,  the  Mg-H2  potential 
was  used,  obtained  by  Jung  and  Gerber  by  SCF/MP-4  ab  initio  electronic  structure  calcu¬ 
lations.  (Ref.  37)  Later,  Chaban  and  Gordon  provided  their  results  for  the  same  systems. 
(Ref.  38)  The  two  potentials,  computed  by  different  methods,  are  quite  similar.  The 
Mg-H2  potential  has  a  slightly  deeper  well-depth  than  the  H2-H2  interaction.  The 
anisotropy  of  the  Mg  —  H2  potential  at  equilibrium  distance  is  very  weak  (weaker  than 
that  for  B  -  H2,  mentioned  in  the  previous  section),  and  can  safely  be  ignored  in  simula¬ 
tions  of  Mg  -  (para)H2. 

(2)  Ground-state  structures  of  Hg(H2)12,  Mg(H2)12:  Figure  7  shows  the  Hg-H2  in  the 
vibrational  ground  state  of  the  latter.  (Ref.  30)  Figure  8  shows  the  H2  -  H2  distance  distri¬ 
bution  for  the  same  system.  (Ref.  30)  Analysis  of  these  results  and  of  other  data  from  the 
DQMC  calculations  shows  that  the  structure  corresponds  to  Hg  being  inside  a  solvation 
shell  of  the  H2  molecules.  (Ref.  30)  The  solvation  shell  is  nearly  spherical.  Figure  9 
shows  the  Mg-H2  distance  distribution  for  the  ground  state  of  Mg(H2)12,  and  Figure  10 
gives  the  H2  —  H2  distance  distribution  for  that  system.  (Ref.  30)  Again,  the  results  show 
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Figure  7 

The  HrHg  Distance  Distribution  in  the  Ground  State  of  HgCH^ 
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Figure  9 

The  H2-Mg  Distance  Distribution  in  the  Ground  State  of  Mg(H2)12 
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Figure  10 

The  H2-H2  Distance  Distribution  in  the  Ground  State  of  Mg(H2)i2  (From  DQMC) 
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Figure  11 

ice  Distribution  of  the  Mg(H2)12  at  T=2K  and  at  T=4K 
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Figure  12 

The  H2-H2  Distance  Distribution  for  Mg(H2)12  at  T=2K  and  at  T=4K 


Scale  of  Figure:  23.8  Angstrom 


Figure  13 

Density  Contours  for  H2  Molecules  in  a  Plane  Containing  the  Center  of  Mg,  for  the  Cluster 

Mg(H2)12  at  T=0K 


the  Mg  atom  to  be  "solvated"  inside  a  "shell"  of  H2  molecules.  (Ref.  30)  There  is  essen¬ 
tially  zero  probability  to  find  the  Mg  at  or  outside  the  solvation  shell,  which  therefore 
shields  very  well  the  Mg  atom.  This  is  already  of  importance  in  assessing  the  stability  of 
Mg  on  solid  hydrogen  against  recombination,  etc. 

(3)  Temperature-dependence  of  the  structure  of  Mg(H2)„:  Figure  11  shows  the  Mg-H2 
distance  distribution  for  Mg(H2)12  at  T=2K  and  T=4K.  Figure  12  gives  the  H2  -  H2  dis¬ 
tance  distribution  in  the  same  calculations.  (Ref.  32)  These  results  were  obtained  by  using 
Feynman  Path  Integral  Simulations.  (Ref.  32)  The  results  show  that  the  structure  of  the 
cluster,  while  somewhat  expanded  in  size,  is  as  in  the  ground  state:  An  Mg  atom  solvated 
inside  a  shell.  The  results  show  that  this  behavior  is  maintained  until  the  evaporation  tem¬ 
perature  (7=4.2K)  of  this  cluster  is  reached.  These  results  suggest  that  in  solid  H2,  Mg 
should  be  "shielded"  by  the  neighboring  hydrogens,  locked  in  its  cage  until  jT=4K  at 
least,  and  negligible  diffusivity  should  be  expected. 

(4)  Anisotropy  of  Mg(H2)„  clusters:  The  structure  of  Hg(H2)i2  is  nearly  isotropic.  Figure 
13  shows  the  density  contour  for  H2  molecules  in  a  plane  that  contains  the  center  of  the 
Mg  atom.  The  result  is  for  T=0K.  Clearly,  the  structure  of  Mg(H2)12  has  an  appreciable 
anisotropy,  with  a  "hole"  in  the  envelope  of  the  H2  molecules.  This  behavior  of  an 
anisotropic  quantum  cluster  is  maintained  until  the  evaporation  temperature  of  the  cluster 
(T=4.2K),  at  least  to  a  large  extent.  For  Mg(H2)16,  Mg(H2)17,  Mg(H2)18,  the  first  solvation 
shell  around  the  Mg  is  complete,  and  the  structure  around  the  Mg  is  closer  to  isotropic, 
although  not  perfectly  so. 

(5)  Vibrational  wavefunction  model,  and  vibrational  properties  of  Hg(H2)12,  Mg(H2)12: 
An  analytical  wavefunction  model  was  fitted  to  the  numerically  calculated  DQMC  distri¬ 
bution.  For  the  n  =  12  clusters,  the  fitted  wavefunction  that  worked  well  was: 

m,  •  •  ■ ,  R12,  R)  =  Jfri.fi  1  <  i*j  z  12/n  Z(n)  (9) 

i=l 

where  R,  is  the  position  vector  of  the  i-th  H2  molecule,  R  -  the  position  vector  of  M, 
rt  j  =  IR,  -  R j  I  is  the  distance  between  the  i-th  and  y-th  molecules,  r,  =  IR,  -  Rl  is  the  dis¬ 
tance  between  H2(z)  and  M.  For  j(r,),  the  vibrational  ground  state  of  the  M  -  H2  dimer 
was  used.  7  is  a  Jastrow  factor  (37)  defined  here  by: 

J{ru\  1  <  i*j  <  12}  =  Cexp[-ZX£V(r,,)]  (10) 

i*j 

Urep(r)  =  (j)n  (11) 

C  is  determined  by  the  normalization  of  the  wavefunction  and  a  is  the  only  free  parame¬ 
ter  used  in  fitting  the  wavefunction  to  the  DQMC  results.  Trial  wavefunctions  of  this  type, 
and  more  elaborate  ones,  were  employed  in  variational  calculations  of  quantum  clusters. 
(Ref.  25)  Note  that  the  procedure  used  here  does  not  involve  a  variational  approximation, 
but  a  fitting  to  DQMC  results.  The  quality  of  the  model,  referred  to  in  this  report  as  "the 
spherical  model"  of  M(H2)12  is  shown  in  Figures  7  and  8  for  M=Hg,  is  excellent.  For 
M=Mg  the  fit  is  less  desirable  than  for  Hg,  but  also  quite  satisfactory.  The  wavefunction 
of  Equations  (10)-(12)  offers  a  simple  physical  picture  for  the  vibrational  dynamics: 
Hg(H2)12  has  two  types  of  vibrations:  "Radial"  vibrations  of  H2  versus  Hg,  and  collective 
vibrations  taking  place  on  the  sphere  of  the  (H2)12  shell.  This  should  be  useful  for  the 
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interpretation  of  future  spectroscopy  on  such  systems. 

The  results  of  this  section  suggest  that  Mg  in  solid  H2,  because  of  the  spontaneous 
solvation  of  Mg  inside  Mg(H2)12,  is  likely  to  prove  a  more  advantageous  material  for 
HEDM  than  Li  in  hydrogen.  It  should  be  possible  to  prepare  the  material  with  higher  con¬ 
centrations  of  the  dopant  in  the  case  of  Mg,  and  the  diffusivity  of  Mg  should  remain  very 
low  under  4K,  preventing  segregation.  These  results  accomplished,  together  with  those 
for  B(H2)n,  Tasks  (1),  (2)  and  (4)  of  the  Program  Plan. 

2.6  Photodynamics  and  Spectroscopy  of  Li(H2)2 

The  purpose  of  this  part  of  the  project  was  to  calculate  the  dynamics  in  time  fol¬ 
lowing  electronic  excitation  of  a  model  of  a  doped  hydrogen  cluster,  and  to  characterize 
its  spectroscopic  consequences.  Electronic  excitation  of  metal  atoms  in  solid  hydrogen 
and  also  in  hydrogen  clusters  has  emerged  as  an  important  experimental  tool  in  the  field, 
and  calculations  for  a  simple  model  system  may  help  guide  future  experiments  with 
regard  to  the  main  effects  anticipated.  The  model  chosen  was  the  collinear  cluster 
Li-H2  -H2.  The  position  of  the  Li  in  the  end,  corresponds  to  the  fact  that  in  real  3D 
Li(H2)„  clusters,  Li  is  outside  the  "droplet"  of  the  H2  molecules.  The  calculations  were 
done  by  a  time-dependent  quantum  wavepacket  method.  Since  the  system  is  a  small  one, 
a  "conventional"  grid  method  was  used,  rather  than  the  novel  CSP  method,  although  later 
tests  have  shown  that  the  CSP-based  algorithm  gives  the  same  results.  The  system  has  X 
and  11-type  excited  states,  and  calculations  for  excitation  into  the  II  state  were  carried 
out.  Below  is  a  summary  of  the  results  from  Reference  (40): 

(1)  There  are  two  predissociation  processes, 

Li*(H2)2  ->Li*-H2  +  H2 
and 

Li*(H2)2  — »  Li*  +  (H2)2 , 

with  the  latter  process  having  a  much  higher  yield. 

(2)  A  long  predissociation  time,  of  the  order  of  f>104  fs,  was  found  for  the  excited  clus¬ 
ter.  Snapshots  in  time  of  the  wavepacket  calculated  for  the  process  are  shown  in  Figure 
14. 

(3)  The  photoproduct  Li*  -  H2  is  formed  with  highly  excited  final  vibrational  state  distri¬ 
butions.  The  result  is  shown  in  Figure  15. 

(4)  A  hybrid  quantum-classical  method  was  tested,  in  which  the  initial  state  is  described 
quantum-mechanically,  but  the  dissociation  dynamics  is  treated  classically.  The  conclu¬ 
sion  was  that  for  some  properties  (see  Figure  15,  for  instance)  this  method  does  very  well, 
but  for  the  lifetimes  and  branching  ratios  it  is  severely  in  error,  and  quantum  calculations 
must  be  used. 

These  results  accomplish  Task  6  of  the  Workplan. 

2.7  Collision  Dynamics  and  Dynamical  Stability  of  Hydrogen  Clusters 

The  survival  probabilities  of  quantum  clusters  are  data  of  importance  for  certain 
potential  approaches  that  have  been  proposed  for  preparing  HEDM  propellants.  To 
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Figure  14c 

Snapshot  of  the  Time  Evolution  of  the  Wave  Packet  of  Li*(H2>2  After  Photoexcitation 

t=7200  fs 


Figure  14d 

Snapshot  of  the  Time  Evolution  of  the  Wave  Packet  of  Li*  (02)2  After  Photoexcitation 

t=14400  fs 
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012345678 
H2-Li  Vibrational  State  (v) 

Figure  15 

Final  Vibrational  State  Distribution  of  the  Dissociation  Fragment  Li*  -  H2  After  the 
Photoexcitation 
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Figure  16 

The  Probability  for  Cluster  Survival  in  Collinear  Collisions  of  Li  Atoms  with 

(H2>2  and  (1)2)2 


Figure  17 

The  Survival  Probability  of  B-H2  in  Collinear  Collisions  with  (H2)  Molecules 
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provide  such  data  quantum  calculations  on  the  collision  dynamics  and  survival  probabil¬ 
ity  in  collisions  for  several  clusters  such  as  (H2)2,  B  -  H2,  (H3)2,  etc.  were  carried  out,  in 
collisions  with  H2  molecules  or  with  Li  atoms.  Figure  16  from  Reference  (15)  gives  the 
survival  probabilities  for  (H2)2)  and  (D2)2  in  collisions  with  Li  atoms,  as  a  function  of  col¬ 
lision  energy.  Surprisingly  high  survival  probabilities  are  found,  even  for  collision  ener¬ 
gies  much  greater  than  the  (very  weak)  driving  energy  of  the  cluster.  Figure  17  gives  sur¬ 
vival  probabilities  for  H2-B  in  collisions  with  H2.  We  used  a  collinear  model,  and  com¬ 
puted  both  the  H2  +  H2  -  B  and  the  H2  +  B  —  H2  types  of  impact.  Again,  very  high  sur¬ 
vival  probabilities  are  found.  Analysis  shows  that  the  high  probabilities  are  due  to  a  quan¬ 
tum  effect:  Configurations  in  the  collisions  that  are  conducive  to  collision  induced  disso¬ 
ciation  are  classically  quite  likely,  but  have  a  much  lower  weight  in  quantum  mechanics. 
Further  analysis  of  the  dynamics  is  given  in  Reference  (15).  The  result  shows  an 
approach  to  making  HEDM  systems  that  depends  on  a  high  probability  of  survival 
(dynamical  stability)  of  quantum  clusters  may  be  viable.  This  result  accomplishes  Task  5 
of  the  (revised)  Program  Plan. 

2.8  Photodynamics  of  Mg  Atoms  in  Solid  Hydrogen 

Although  this  part  of  the  work  is  not  completed  yet,  the  results  already  obtained  on 
this  topic  are  important,  and  contribute  to  the  success  of  the  project.  Preliminary  calcula¬ 
tions  on  the  S  — >  P  excitation  of  magnesium  atoms  in  solid  hydrogen  were  performed. 
The  calculations  were  carried  out  by  using  the  time-dependent  quantum-mechanical 
method,  based  on  the  CSP  approximation,  that  was  developed  in  the  present  project.  A 
simulation  box  that  included  54  hydrogen  molecules  with  periodic  boundary  conditions, 
at  the  density  of  hydrogen  at  normal  pressure  for  T  — » OK  was  used.  The  system  was 
"doped"  with  Mg  atoms  by  removing  4  hydrogen  molecules  from  the  center  of  the  box, 
and  inserting  an  Mg  atom  instead.  Estimates  indicate  that  the  replacement  of  4  hydrogen 
molecules  by  an  Mg  atom  is  energetically  the  most  favorable.  Diffusion  Quantum  Monte 
Carlo  (DQMC)  on  this  system  of  one  Mg  surrounded  by  50  hydrogens,  and  with  periodic 
boundary  conditions  were  then  carried  out.  A  model  wavefunction  was  then  fitted  to  the 
DQMC  results,  along  the  lines  of  the  model  wavefunction  fitting  for  Mg(H2)12,  described 
in  section  5.  To  model  the  S  — »  P  excitation,  the  above  ground  state  wavefunctions  were 
used  as  the  initial  state,  and  began  to  propagate  it  on  the  excited  potential  surface.  Excita¬ 
tion  into  the  singlet  state  was  considered.  Tests  have  shown  that  the  simulation  box  with 
50  hydrogen  molecules  is  too  small,  and  the  results  can  therefore  not  be  trusted,  at  least 
quantitatively.  There  are  indications  that  for  a  simulation  box  with  100  hydrogen 
molecules  the  size  of  the  box  will  be  sufficient  for  physically  realistic  results,  at  least  for 
short  times.  The  qualitative,  or  semiquantitative,  results  that  nevertheless  merit  attention 
and  may  remain  valid  for  the  larger  simulation  box  are: 

(1)  Following  excitation,  photoinduced  mobility  of  the  Mg  is  found  with  high  yield: 
Nearly  35%  of  the  Mg*  atoms  migrate  from  the  original  cage  to  a  neighboring  cage. 

(2)  For  the  Mg*  atoms  that  have  not  migrated,  highly  excited  cage  vibrations  take  place. 
These  vibrations  are  mostly  "radial"  with  the  H2  molecules  oscillating  against  the  Mg. 
The  vibrations  are  coherent  for  r~400  fs  at  least,  and  decay  slowly  due  to  coupling  with 
other  H2  molecules. 

(3)  Appreciable  non-adiabatic  transitions  are  seen  already  for  t>  600  fs.  In  these  events, 
the  electronic  population,  initially  in  the  highest  P-state,  relaxes  in  part  to  the 
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intermediate  electronic  state.  There  is  no  significant  electronic  population  in  the  lowest 
electronic  state,  at  least  for  t<  800  fs. 

It  should  be  extremely  useful  to  study  the  system  of  Mg  in  solid  hydrogen  by  ultra¬ 
fast  time-resolved  spectroscopy,  to  observe  both  the  photo-diffusion  and  the  nonadiabatic 
transitions.  Data  on  electronic  spectroscopy  of  Mg  in  solid  H2  was  reported  by  Presser. 
(Ref.  41) 

These  results  contribute  to  achieving  Task  2  of  the  Program  Plan  of  our  research, 
and  with  them  also  Task  8  is,  in  part,  accomplished. 

2.9  Structural  and  Thermodynamic  Properties  of  Hydrogeh/Rate  Gas  Mixtures  on 
Surfaces 

This  result  touches  only  in  a  marginal  way  on  our  project  and  also  there  is  only  a 
very  preliminary  result.  Nevertheless,  because  of  the  scientific  interest,  it  is  mentioned 
here.  It  was  suggested  almost  three  years  ago  that  monolayers  which  are  mixtures  of  H2 
with  Rg,  should  have  disordered  structures  that  are  very  different  from  corresponding  rate 
gas  mixtures.  For  instance,  mixtures  of  Xe+Kr  on  Pt  are  substitutionally  disordered,  but 
crystalline  -  periodic  over  a  wide  range  of  T.  For  H2/Ar  mixtures  on  a  flat  surface  we 
expected  other  types  of  disorder,  and  a  strong  effect  of  T  on  the  disordered  type.  This  was 
due  to  the  fact  that  for  Rg,  the  van  der  Waals  radii  play  the  critical  role  and  determine  the 
disordered  type,  while  for  H2,  the  de  Broglie  wavelength  is  the  key  parameter. 

A  cooperation  was  established  with  P.  Zeppenfeld  (Julich))  and  M.  Bienfait  (Mar¬ 
seilles)  which  carried  out  neutron  diffraction  experiments  on  H2+Ar  on  graphite,  in 
Grenoble,  France.  Early  results  suggest  indeed  a  possibility  of  a  novel  type  disorder 
behavior. 

The  results  may  be  useful  in  understanding  the  structures  and  thermodynamic  sta¬ 
bility  of  H2+M  mixtures  (M-doping  material)  and  high  mixing  ratios. 

3.  WORK  IN  PROGRESS 

Calculations  are  in  progress  on  the  dynamics  following  electronic  excitation  of  Mg 
atoms  in  solid  hydrogen.  Results  were  reported  in  Section  (2.8),  but  the  number  of  H2 
molecules  used  in  the  simulations  was,  by  the  current  evidence,  too  small.  Simulations  of 
Mg  in  a  sample  of  over  108  hydrogen  molecules  (with  periodic  boundary  conditions) 
should  be  necessary  for  realistic  physical  results.  The  calculations  carried  out  now  are 
time-dependent  quantum  mechanical  simulations,  by  the  method  developed  in  this  pro¬ 
ject.  Properties  such  as  photoinduced  mobility,  electronic  energy,  vibrational  relaxation 
(of  the  excited  "spot"  consisting  of  the  Mg  atom  and  the  surrounding  cage),  and  elec¬ 
tronic  adsorption  spectroscopy  are  being  computed. 

Another  calculation  still  in  progress  is  the  continuation  of  the  calculations  for 
B(H2)„  clusters,  with  n  -  1, . . . ,  8.  Currently,  DQMC  calculations  of  B(H2)„  with 
n  ~  100  are  being  performed  along  with  calculations  for  B  in  solid  H2,  with  periodic 
boundary  conditions.  The  aims  of  the  calculations  are  structural  properties,  energetic  sta¬ 
bility,  and  consequences  for  diffusion.  Calculations  of  electronic  spectra  are  also  being 
contemplated. 
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4.  RECOMMENDATIONS  AND  CONCLUSIONS  FOR  THE  HEDM  PROJECT 


Included  here  are  only  a  few  main  recommendations.  More  detailed,  technical 
conclusions  and  recommendations  are  mentioned  in  the  Summary  and  in  the  Results 
Sections. 

(1)  Focussing  with  high  priority  on  developing  boron  in  solid  hydrogen  as  a  HEDM 
propellant.  By  the  aspects  studied,  this  is  by  far  the  most  promising  HEDM  candidate. 

(2)  Mg  in  hydrogen  seems,  by  the  results  of  this  research,  to  have  important  advantages 
over  Li  in  hydrogen  as  a  potential  HEDM  propellant.  It  appears  that  now  Li  in  H2  is 
getting  more  attention.  By  the  reported  results,  Mg  in  H2  has  better  stability  and  can  be 
built  up  to  higher  doping  concentrations.  This  suggests  a  change  in  priority.  Boron  in 
hydrogen  seems,  nevertheless,  much  more  promising  than  either  of  these. 

(3)  The  Phillips  Laboratory  should  consider  time-resolved  spectroscopy  as  a  new, 
additional  tool  for  in-house  experiments  to  characterize  HEDM  systems,  especially  as  to 
dynamical  properties,  and 

(4)  The  Phillips  Laboratory  should  consider  adopting  time-dependent  quantum 
simulations  as  a  tool  for  supporting  some  of  the  experiments. 
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